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Abstract. We study the ground-state properties of ultracold bosons in an optical 
lattice in the regime of strong interactions. The system is described by a non- 
standard Bose-Hubbard model with both occupation-dependent tunneling and on-site 
interaction. We find that for sufficiently strong coupling the system features a phase- 
transition from a Mott insulator with one particle per site to a superfluid of spatially 
extended particle pairs living on top of the Mott background - instead of the usual 
transition to a superfluid of single particles/holes. Increasing the interaction further, a 
superfluid of particle pairs localized on a single site (rather than being extended) on top 
of the Mott background appears. This happens at the same interaction strength where 
the Mott-insulator phase with 2 particles per site is destroyed completely by particle- 
hole fluctuations for arbitrarily small tunneling. In another regime, characterized by 
weak interaction, but high occupation numbers, we observe a dynamical instability in 
the superfluid excitation spectrum. The new ground state is a superfluid, forming a 
2D slab, localized along one spatial direction that is spontaneously chosen. 
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1. Introduction 

Systems of ultracold atoms in optical lattices provide a unique playground for controlled 
realizations of many-body physics [TJ [2]. For sufficiently deep lattices, the kinetics 
is exhausted by tunneling processes, and an initially weak interparticle interaction 
eventually becomes important with respect to the kinetics, when the lattice is ramped 
up. A consequence of this competition is the quantum phase transition from a superfluid 
of delocalized bosons to a Mott insulator where the particles are localized at minima 
of the lattice by a repulsive contact interaction [3]. This effect has been observed in 
seminal experiments with ultracold rubidium atoms in a cubic lattice pE]. It is described 
quantitatively by means of the simple Bose-Hubbard model [3], [5], whose parameters are 
the interaction energy U for each pair of particles occupying the same lattice site, and 
the matrix element J for tunneling between neighboring sites. Intriguing Hubbard-type 
physics can also be observed if the above scenario is extended to fermions, mixtures 
of several particle species, exotic lattice geometries, or long-ranged dipolar interaction 

DO 121 IB] . 

In this paper, we consider a different type of extension of the bosonic Hubbard 
model, becoming relevant when the interaction between the particles is enhanced, e.g, 
by means of a Feshbach resonance. As long as the interaction is weak compared to the 
lattice potential, a system of ultracold atoms can be described, to a good approximation, 
in terms of the lowest-band single-particle Bloch or Wannier states, the latter being 
localized at the minima of the lattice [7]. Under these conditions, the Hubbard 
interaction U and tunneling parameter J are given by respective matrix elements 
with respect to the single-particle Wannier states. This approximation corresponds 
to degenerate perturbation theory up to first order with respect to the interaction, 
with only the intraband coupling induced by the interaction taken into account. 
However, if the interaction is stronger, higher-order corrections start playing a role. 
One may still describe the system in terms of lattice-site occupation numbers rij, but 
the occupied Wannier-like orbitals will have admixtures from higher bands, depending 
on the occupation. The most significant effect of the repulsive interaction will be a 
broadening of the Wannier-like orbitals with increasing occupation, effectively enhancing 
J and decreasing U. In terms of the Hubbard description, we take this into account by 
replacing J and U by functions Jh^hj and Uf H of the number operators hi. Quantitative 
consequences of this kind of modification to the plain bosonic Hubbard model have been 
studied by several authors at a theoretical level [H HI [10]. Considering an interaction- 
induced modification of the Wannier functions, also additional Mott-insulator phases 
have been predicted [TT1 IT2] . In Ref . [13] , the effect of the interaction- induced coupling 
to the first excited band on the Mott transition was considered. Re-entrant behaviour 
in the superfluid-Mott transition has also been predicted due to the interaction-induced 
modification of Hubbard parameters [131 E3]- The effect of interaction on the tunneling 
dynamics in one-dimensional double-well and triple-well potentials have been studied 
in Refs. [151 EE] where the authors found enhanced correlated pair tunneling near the 
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fermionization limit. Moreover, occupation-number-dependent on-site interaction has 
been observed experimentally in the coherent dynamics of an atomic ensemble |17j . 
Similar occupation-dependent effects have been observed in Bose-Bose [18] and Fermi- 
Bose mixtures [19] [20], [21], and - in the latter case - have been explained theoretically 
in terms of occupation-dependent parameters U and J [22] . 

In the present work, we show that new quantum phases can arise in Hubbard 
models with number- dependent parameters. After writing down the effective single-band 
Hamiltonian including the effect of the site occupation, we find that for strong enough 
interaction (characterized by the s-wave scattering length a s ), there is a transition from 
a Mott state with one particle localized at each lattice site to a superfluid of pairs 
extended over neighboring sites, rather than to a superfluid of single atoms. This 
feature is novel, considering the fact that the extended pairs emerge in the single- 
species repulsive bosonic system without the presence of any long-range interaction. 
For even higher interaction strengths, the n = 1 Mott state becomes unstable towards 
a superfluid of pairs that are localized on single sites. Moreover, the n = 2 Mott state 
becomes unstable towards pair fluctuations already for very low tunneling amplitudes. 
Finally, we consider the regime where interaction effects are important not because of 
large scattering lengths, but rather because of large site-occupations numbers. In this 
limit, starting from the Bogoliubov approach to the homogeneous system, we find a 
phonon instability at a critical filling fraction. Above that fraction, the new ground 
state is a Bose condensate with the particle density being localized along one spatial 
direction that is chosen spontaneously. 

This paper is arranged in the following way: In section [2] we introduce the 
occupation-dependent Bose-Hubbard model. In section [3] we start discussing the 
properties of this model. Namely, we study the instability of the Mott-insulator phase 
with respect to simple particle and hole excitations, leading to the usual single-particle 
superfluidity. In section SJ we then investigate the instability of the Mott phase 
with respect to the excitation of bond-centered pairs of particles being extended over 
neighboring lattice sites. We show that this mechanism will eventually become relevant 
when the s-wave scattering length is increased, and that one finds a phase transition 
to a superfluid of extended pairs. In section EJ proceeding to even stronger interaction, 
the instability of the Mott phase towards a superfluid of site-centered pairs is discussed. 
In this regime, moreover, the Mott insulator at a filling of two particles per site can 
disappear completely. Finally, in section EJ we focus on the limit where interaction- 
induced orbital effects play an important role because of large filling. We find that, with 
the increasing superfluid density, the condensate may become dynamically unstable. 
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2. The Bose-Hubbard model 



The Hamiltonian in the presence of a periodic potential with lattice constant a, given 



by K 



per V 



Vo[sin 2 (7re/a) + sm 2 (ny/a) + sin 2 (irz/ a)], reads 
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with bosonic field operators ^, mass m, and interaction strength g = A7rh 2 a s /m, where 
a s is the s-wave scattering length. To derive a Hubbard-type description, the field 
operators V>(r) are expanded in terms of Wannier-like orbitals <f>i(r,hi) = 4>(r — Ri,hi) 
localized at the lattice minima Ri, namely ifj(r) = ^2ibi<f)(r — R^hi) with bosonic 
annihilation and number operators hi and hi = b\bi. Note that the "wave function" (pi 
depends on the number operator hi in order to take into account interaction-induced 
occupation-dependent broadening. Keeping only on-site interaction, we arrive at the 
effective single-band Hamiltonian 
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and we have introduced the chemical potential \i to control the particle number. We 
would like to mention that in the presence of an optical lattice for high interactions the 
pseudo-potential form of contact interaction can still be used, when a modified scattering 
length which is different from the bare scattering length is applied (23J EU [251 126] . 

In order to estimate the occupation number dependence in a mean-field way, 
we make a Gaussian ansatz for the Wannier-like wave functions, <f)(r — Rf, n^) = 
exp(— (r — R) 2 /d 2 (ni)), where the width d(ni) is a variational parameter depending on 
the particle number rij, and minimize the Gross-Pitaevskii energy functional. The idea 
to use the width of the Wannier function as a variational parameter has also been used 
in Refs. [271 12H1 [29] . Taking into account the full lattice potential (i.e., not employing 
a quadratic approximation for the lattice minima), for a given rij this leads to 
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/it for the width of in the limit Vq ^> Er, where 



Er = w 2 h 2 /2ma 2 denotes the recoil energy. Note that Eq. (j4j) has a solution only as 
long as \/Vol ' Er ^> d 2 {ni)/d^. Using the variational result, the tunneling parameter 
between two adjacent sites can be approximated by 
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We would like to point out that when calculating the tunneling strength, the Gaussian 
approximation generally results in a lower value than the exact calculation; the 
exact Wannier orbital has an exponential tail which decays slower than a Gaussian. 
Nevertheless, our simple approximation provides us with reasonable numerical values 
and with a suitable model for the occupation dependence of tunneling in the regime 
treated here. This allows us to get a qualitative understanding of the physics at work. 
For the number-dependent on-site interaction strength the variational result gives 



3/2 

do 



Er v "\E r J ] L {d 2 {n i ) + d 2 {n i -l))\ a' (6) 
The single-particle tunneling term arising from the non-on-site contributions of 
the quartic interaction term in Eq. ([T|) is exponentially smaller than J(ni,Uj) by 
approximately a factor of exp(— it 2 ^Vq/ ErJ A)a s / a. Similarly the pair tunneling term 
is smaller than J(ni,nj) by approximately a factor of exp(— n 2 ^Vq/ Eji/2)a s /a. Since 
we are in the limit of Vo/Er 3> 1, these terms are neglected in Eq. ((2]). 

3. Insulator to single-particle superfluid transition 

Having written down a suitable model Hamiltonian describing the regime of strong 
interaction, we now study the transition from the Mott insulator having on average n 
particles per site to a superfluid of single particles/holes. 

For this purpose, we use a product ansatz Y\i f° r the many-body state, with 
the variational coherent spin-representation state [301 l3Tj . 

|<E>)i = cos#|n)j + sin 9 sin^|n + 1); + sin 9 cosip\n — l)j (7) 

at each site i, with occupation number basis states In*)*. Here we only take into 
account states with one additional particle or hole, which in the Mott phase and close 
to the transition to the superfluid, where particle fluctuations are small, is sufficient. 
Accordingly, the variational mean-field energy is given by 
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where 



Hj = (n 2 + n)J % n sin 2^/2 

+ (n + 1)7^+1^ sin 2 ?/; + nJ^n-i cos 2 ip, 
Hu = n(n — l)L^:Cos 2 9 + nin + l)Un+i sin 2 #sin 2 tp 

+ (n- l)(n-2)f4„ 1 sin 2 0cos 2 ^. (9) 

Minimizing the energy determines 9 and ip. While 9 = corresponds to an 
incompressible Mott-insulator state with an integer number of particles n per site (found 
within a finite interval of the chemical potential /x), the superfluid state is characterized 
by 9 7^ with order parameter ~ sin 29. In the superfluid phase, the average particle 
number per site is characterized by if) depending smoothly on the chemical potential. For 
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ip <C 7r/4, the transition to the superfluid occurs mainly via the creation of holes, while 
for ip near ir/2 particle creation is the main mechanism destroying the Mott phase. In 
the latter case, the Mott insulator becomes unstable when the energy cost of creating an 
additional particle at one site, namely U n+ in(n + 1) / 2 — /i, is overcome by the reduction 
in energy due to tunneling of that particle, which is on the order of z(n + 1) Jn+i,n, with 
coordination number z = 6 for the cubic lattice. Thus, when E ss minimizes for non-zero 
6, the Mott state becomes unstable with respect to single particle and hole excitations. 
For interaction strength a s /a = 0.15 and n = 1 , this happens at the black lines (solid 
or dotted) in the plane spanned by /x/Vq and Jo,i/V in Fig. [TJ 




Figure 1. Mott-insulator-to-supcrfluid phase transition for a s /a = 0.15. Inside the 
region marked by the black solid line and the blue dashed line, the system is a Mott- 
insulator with n = 1 particles per site. Leaving this region by crossing the black solid 
line, a simple superfluid of single particles (or, equivalently, holes) is formed (SF). In 
contrast, crossing the blue dashed line one arrives at a superfluid phase of extended 
(bond-centered) pairs (ePSF). In technical terms of our variational approaches: outside 
the black solid and dotted line minimizing the energy ([8]) gives 6* =^ 0, while on the 
r.h.s. of the blue dashed line 9 e ^ is obtained from minimizing expression (|12[) . 



4. Superfluidity of extended (bond-centered) pairs 

So far, we have described the usual scenario of the Mott phase becoming instable with 
respect to particle and hole derealization, as it is also found for non-number-dependent 
Hubbard coupling J and U. However, we will now show that — as a consequence of 
occupation dependent hopping and on-site interaction — the Mott insulator with n = 1 
can become unstable with respect to the creation of pairs of particles already before 
the creation of single particles becomes favorable. Consider a pair excitation with one 
additional particle at site i and another one at the neighboring site j, corresponding to 
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the state \P(ij)) = \b\bj\{ n i — !})• Such a bond-centered or extended pair excitation at 
(ij) can tunnel coherently to a neighboring bond, say (ik), with k ^ j being another 
neighbor of i. Generally, bonds are considered neighbors if they share a common 
site. Such a pair tunneling processes occurs in second order with respect to single- 
particle tunneling via the virtual site-centered pair state |Pj) = -^b\b\\{rii = 1}}, 
which has larger energy. According to second-order degenerate perturbation theory, 
the amplitude of the pair tunneling process is given by J c ^ = 6J| 2 /(3f/ 3 — 2U 2 ). On 
the same footing, perturbation theory gives the binding energy of —2J e g for the bond- 
centered pair due to number fluctuations within the pair. For a cubic lattice of sites, the 
bond-centered pair excitations live on an exotic lattice of coordination number z' = 10, 
being a generalization of the two-dimensional checkerboard lattice (see the rightmost 
drawing in Fig. [2]) to three dimensions. This allows the pair to reduce its energy by 
10J c fj when delocalizing. In contrast, two additional particles, not forming a pair, can 
reduce their energy by 2 x 6 x 2J\ 2 when delocalizing on the cubic lattice of sites 
(coordination number 6). Thus, according to perturbation theory, the formation of a 
bond-centered pair is favorable if —(10+2) J c g > 24J 12 . For certain scattering lengths a s , 
this condition can be fulfilled, since the Wannier-broadening with increasing scattering 
lengths leads to an increase of both J 2j 2/^1,2 an d U 2 /U 3 . In such a situation, the Mott- 
insulator state becomes unstable with respect to the creation of bond-centered pairs 
rather than with respect to the creation of single-particle excitations. This happens 
when the derealization energy — 10J c fj overcomes the energy 2(U 2 — /1) — 2J cS needed 
to create a pair excitation. It is interesting to note that an equivalent scenario does not 
happen for hole excitations, since hole excitations decrease the occupation number and 
with that the tunneling amplitudes. 




Figure 2. Color online: The left hand side shows a square lattice of sites (blue squares) 
connected by bonds (black lines). The lattice of the bonds of the square lattice, where 
bonds sharing a site are connected, is given by the checkerboard lattice shown on the 
right hand side. If a bound pair of two indistinguishable particles can either occupy 
a site or a bond of the cubic lattice (the latter means that the two particles occupy 
neighboring sites) and if the pair can move (by single-particle tunneling) from a site 
to a neighboring bond and vice versa, then the pairs move on the lattice shown in the 
center plot. Sites and bonds arc denoted by blue squares and black bullets, respectively. 
Extending all the considerations shown in this figure to the case of a three-dimensional 
cubic lattice of sites is straightforward. 
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To evaluate the boundary of the n = 1 Mott-insulator phase within mean-field 
approximation, we construct a model for the excited bond-centered pair excitations. 
When the number of pairs is small compared to the number of sites, the Hamiltonian 
for the pairs living on top of a Mott state with one particle per site can be written as 



Here, L = labels the bonds of the cubic lattice and (LU) denotes pairs of nearest 
neighbors of these bonds as they are described by the three-dimensional checkerboard 
lattice (cf. Fig. |2]). Moreover, we have defined the bosonic creation and destruction 
operators for bond-centered pair-excitations p' L and pi, with number operator hi = 
P l Pl- As a consequence of the diluteness assumption, we have neglected the interaction 
between pairs, arising if pairs occupy neighboring bonds. Since the transition to a pair- 
superffuid will happen with the creation of a single pair, this approximation will not 
influence the phase boundary. The energy of a condensate of bond-centered pairs can 
now be estimated in a similar fashion as before by making a product ansatz FJ L \& p )l 
of coherent states being a superposition of zero and one pair at each bond, 



where z' = 10 is the coordination number of the three-dimensional checkerboard 
lattice. The mean-field approach gives the same phase boundary for the appearance 
of a pair condensate with finite order parameter (pi) as the perturbation theoretical 
considerations of the previous paragraph. The equivalence of both approaches is 
generally given for an ansatz like OTp which includes only two states per site. 

In Figure [TJ we plot the results of minimizing E ss , E ep with respect to 9, 9 e for 
a s /a = 0.15. The stable Mott region with respect to single particle-hole excitation 
is given by the interior of the black solid and dotted line characterized by 9 = 0. 
On the right hand side of the blue dashed line in Fig. [TJ one finds a region where 
min[_E cp ] < min[_E ss ] with 9 e ^ 0. Thus, here the system is characterized by (p£) ^ 
and (fej) = 0, i.e., the state is a superfluid of extended pairs (ePSF). 

Condensates of extended pairs have also been proposed in the context of dimer 
models of reduced dimensions, describing frustrated magnets like SrCu2 (603)2 [32]. By 
approximating triplet excitations as hard-core bosons, the authors of Ref. [33] argue that 
for correlated hopping these bosons can condense in pairs. Such pairing processes also 
bear resemblance to molecular condensation due to Feshbach resonances in an optical 
lattice [33| . 

We would like to point out that triple, quadruple or higher order excitations do 
not play a dominant role. The effective tunneling matrix element of such excitations 
will be very small since it appears in third or higher order perturbation theory only. 
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Therefore inside a n = 1 phase, triple and higher excitations cannot lower their energy 
efficiently by derealization. We can, thus, exclude a superfluid of triples or higher order 
objects. However, there is another possible and competitive scenario we would like to 
mention. Instead of exciting a triple or quadruple, one can create a huge cluster of extra 
particles, i.e., a big spatial domain with doubly occupied sites. In this case, within each 
cluster, the energy of the additional particles (on top of the n = 1 Mott background) 
is not lowered by derealization, but rather by the attractive interaction between them 
as it appears in second order perturbation theory. In the bulk of such a cluster, this 
gives a binding energy of —QJ e s per extra particle. In comparison, in the pair superfluid 
each particle can lower its energy by J c s because of binding and further by another 
5 J e ff because of derealization (i.e., Bose condensation). Accordingly, in leading order a 
superfluid of bond-centered pairs on top of the n = 1 Mott insulator is equally favorable 
as a phase separated state with spatial domains hosting either a Mott insulator of filling 
n = 1 or n = 2. As a consequence, we cannot reliably exclude phase separation by 
means of simple variational arguments. 

Before moving on, let us briefly discuss another issue: In this article, we are working 
in a situation with the chemical potential fixed rather than the particle number. This 
approach is actually quite suitable for the description of experiments with ultracold 
atoms, provided the atoms are trapped by a sufficiently shallow potential. In such 
a situation, the local density approximation applies and different regions in the trap 
correspond to different values of the chemical potential. However, if the trap is too 
steep for the local density approximation to be valid, it might introduce also new 
physics. Consider the following example: The phase separated state described in the 
preceding paragraph might not be favored in the homogeneous system. But, because 
it is energetically very close to the pair superfluid, it can be favored already when a 
slight potential difference is introduced, helping to form n = 2 Mott domains in the 
region of slightly lower potential energy. Such a scenario can spoil the local density 
approximation already for a very weak trapping potential. 

5. Superfluidity of local (site-centered) pairs 

Now, let us consider a regime that can be achieved if the interaction strength a s /a is 
increased further. Considering again the n = 1 Mott insulator, for increasing interaction 
a site-centered pair excitation, described by |Pj), eventually becomes more favorable 
than the bond-centered excitations described by \P(ij))- This happens when the ratio 
U 3 /U2 is reduced so much that the potential energy 3U 3 needed to create a pair of 
particles on the same site equals the potential energy 2U2 required to create a pair 
of particles on neighboring sites. Such a situation is possible as can be derived from 
Eq. (J4]). In the limit of large Vq Er and a s /a we can write d(n)/do ~ (fi 1 ^) 1 ^ 5 
resulting in 3U3 — IXJ2 ~ — 0.02f/o. If \iU 3 — 2C/2 1 becomes comparable to or smaller 
than J 2 2 , a bond centered pair excitation \P(ij\) can transform to a site-centered pair 
excitation |Pj) by a single-particle tunneling process described by the matrix element 
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-Tpair = \f§J 22 - in this regime, the pairs occupy the lattice given by both the sites and 
the bonds of the cubic lattice (see Fig. [21 center). By delocalizing on this lattice, a pair 
can reduce its kinetic energy by 12 J pa ir- As long as this energy is bigger than the kinetic 
energy reduction 24J 12 which two non-paired particles can achieve by derealization, 
the pair is stable towards breaking; this is the case for J 2j2 > \/3/2Ji,2- Thus, the 
binding mechanism of the pair is based solely on the derealization of its center of mass. 
At \3U 3 - 2U 2 \ » (given, e.g., for a s /a « 0.21 when V /E R » 16), the n = 1 Mott 
insulator becomes unstable with respect to pair creation when 12 J pair exceeds 3U 3 — 2/j,. 
It is fascinating to observe the emergence of exotic lattice geometries as illustrated in 
Fig. [2] as a consequence of pair creation. 

If the scattering length is increased further, such that 2U 2 — 3U 3 ^> J 2j2 , site- 
centered pair excitations \Pi) will be created rather than bond centered ones \P(ij))- 
The site-centered pair excitations can then tunnel from site to site coherently via the 
occupation of a virtual bond-centered pair excitation. The corresponding tunneling 

j2 

matrix element reads J' eS = 6 2 u 2 -3U 3 = ~ ^ cS ' Moreover, the pair has a binding energy of 
6 J' cS (stemming from a small perturbative admixture of the 6 neighboring bond-centered 
pair states). Therefore, a site-centered pair is more favorable than two single-particle 
excitations if 3Uz — 12J' eS < 2{U 2 — 12Ji t 2)- If this condition is fulfilled, the Mott insulator 
becomes rather unstable towards the creation of site-centered pair excitations than to 
the creation of single particles. The instability occurs when \2J' cS reaches 3U 3 — 2\i. 
As before, a mean-field calculation leads to the same phase boundary. We plot the 
boundary of the n = 1 Mott phase for a s /a = 0.3 in Fig. EJ The instability towards 
the creation of single particles is hardly important. It is predominantly the creation of 
single holes or site-centered pairs of particles which destroys the Mott phase. 

Note that in the limit of U\ ^> Jq,i metastable repulsively bound pairs of ultracold 
bosons have been observed in optical lattices [351 EE]- Also, two-species mixtures of 
bosons with inter-species attraction trapped in an optical lattice have been shown to 
give rise to superfluidity of pairs [3?]. In the context of dipolar atoms in a two- leg 
ladder, when no tunneling is present between the two legs, pair superfluidity arises due 
to attraction between the dipolar atoms between the two legs of the ladder [381 EH]- 
Also, using a state-dependent optical lattice potential, it is possible to create correlated 
tunneling of on-site pairs, which in turn gives rise to superfluidity of local pairs jlQl E] . 
In our present study, we find that such local pairing can emerge due to the strong 
occupation-dependence of tunneling and on-site interaction. 

After having studied the boundaries of the Mott-insulator phase with one particle 
per site, let us have a look at the n = 2 Mott state. In the limit of vanishing tunneling, 
a Mott state with two particles localized at each site is favorable for U 2 < fi < 3U 3 — U 2 . 
The upper border of this interval is given by the potential energy difference of having 
three and two particles at a site. This difference can, in fact, become lower than 
the potential energy difference U 2 between two and one particle per site marking the 
lower border. This is the case if 3U 3 — 2U 2 < 0; then the n = 2 Mott-insulator 
phase is never stable with respect to the creation of particle-hole pairs, irrespective 
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Figure 3. Mott-insulator-to-superfluid phase transition for a s /a = 0.3. Inside the 
region enclosed by the black solid and the blue dashed line, the system is a Mott- 
insulator with n = 1 particles per site. Crossing the dashed blue line, one enters a 
superfluid of local, site-centered pairs (PSF). Leaving the Mott-phase by crossing the 
black solid line a superfluid of single particles (or, equivalently, holes) is found. Black 
dashed line defined as in Fig. [TJ 



of the tunneling strength; it ceases to exist. The disappearance of the n = 2 Mott 
insulator coincides with site-centered pair excitations becoming more favorable than 
bond-centered ones in the limit of vanishing tunneling. Note that the Mott-insulator 
phases with higher filling, n > 3, do not disappear for large interaction a s /a within 
the Gaussian approximation. The reason why these phases do not share the fate of the 
n = 2 Mott insulator is that the broadening on the Wannier-like site-wave functions <pi 
in response to adding one particle to that site becomes less pronounced with increasing 
occupation: U2/U3 > U3/U4 > > • • ■. However, one should have in mind that for 

strong interaction, sites occupied by three and more particles suffer strong dissipation 
due to three-body collisions [121 1^3] . 

One might ask about the nature of the system's ground state at fixed filling n = 2 
and for 3U3 — 2U2 < 0, when there is no n = 2 Mott phase. At vanishing tunneling, the 
ground state is highly degenerate consisting of all Fock-states having occupation n t = 1 
on half of the sites and occupation rij = 3 on the others. Alternatively, one might say 
that on top of a n = 1 Mott insulator, half of the sites are occupied by additional 
site-centered pairs. For small but finite hopping this degeneracy will be lifted. One 
can think of three possible scenarios: (i) The pairs gather in one region in space; this 
corresponds to a phase segregation between the n = 1 and the n = 3 Mott phases, 
(ii) The pairs delocalize to form a superfluid. (iii) The pairs form a checkerboard-type 
insulator avoiding pairs on neighboring sites. In order to decide this question, we write 
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down an effective Hamiltonian for the site-centered pairs: 

tfpair = -^eff E ^ ~ E^ ~ 6 J '<*K + ( J 'eS ~ A ) E ^ ( 13 ) 

(».J> » (*i> 

with bosonic pair annihilation and creation operators q, ct, and where we assume a 
hard-core constraint (cj) 2 = 0. The nearest-neighbor repulsion present in the last term, 

j2 

with A = 2 6U4+ U2_ 6U3 , stems from super-exchange processes between neighboring pairs. 
This model can be mapped to a spin-1/2 XXZ model with the first term corresponding 
to the XX coupling and the last one to the Z-coupling. Since ( J' cS — A) < J' e{{ is always 
true, the system will neither form the checkerboard pattern (iii) (corresponding to an 
antiferromagnetic state for the XXZ-magnet) nor show phase segregation (i) [33]. The 
system forms a superfluid of site centered pairs (ii). 



6. Weakly interacting limit 



Finally, we investigate the limit where interaction effects are important not because of 
a large scattering length but because of large site occupation, i.e., a s /a <C 1, but the 
mean number of particles per site uq ^> 1. We assume small on-site number fluctuations 
Sn <C n , i.e., ^/U no / (n J no ) 1. In this limit, we can write the modified Hubbard 
Hamiltonian as 



H = - J no E $[1 + a(5hi + Sh^bj 



+ 



E 



nArii 



l)[l + /3-2/3(n i -l)]-E^ 



where 



P 



a 



J, 



3 hi 
5V2 

7T 5 / 2 



e r 
v 



V 



10x/2 [E R 

7T 2 _ 

4 



1/4 



3/4 



a 
a 



cxp 



Er 



Vn 



nl/4 



R 



a 



-n 



(14) 

(15) 
(16) 
(17) 



and 5n 



Hi 



Uq. Here, we would like to point out the similarity of Hamiltonian (il4"|) 

5], given by 



to the Quantum Ablowitz-Ladik (AL) model for g-deformed bosons 

- QBjBi)}, 



E[^+i + 4+i^ + ^ln(l 



where [Bi, Bj] = exp[— 2'jNi], and Q — 1 — exp[— 27]. In the limit of 7 — > and 
7A r j 1, Eq. ( ITS]) reduces to the occupation-dependent modified Hubbard model 
Eq. (|T4|) with a = 7 and U m = 0. It is found that in one and higher dimension the AL 
model contains localized solutions [4*61 Wf\ . To investigate this possibility, we first solve 
Eq. ( !T4|) in the superfluid limit, where the order parameter reads = y/no. To look 
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for fluctuations around the ground state, we first convert the Hamiltonian in Eq. (jl4p 
to momentum space by defining bj = ^ b k exp(— ik.fi), £k = 4^j =12 3 sm2 ^) an d 
7fc = 4 5^i=i 2 3 cos2 ^T- Neglecting correlations arising from the three-body interaction 
term in Eq. ( fl4|) . one arrives at the Hamiltonian 

# m od = " ^ + J n^kb\b k (19) 



2 

k 



^=<L (1 + 0- 2(3(n - 1)) - aJ no n olk 



x (2b\b k + b\b ] _ k + 6 fe 6_ fe ) . (20) 

It can be diagonalized via a Bogoliubov transformation, and the excitation spectrum 
Qk of the superfluid is found to be given by 

Ql = J 2 no e 2 k + 2U no n (l + p- 2f3(n - 1) - 2a^ 7 fc) e k . (21) 

In a cubic lattice, as k — > one finds flk/J no Uo = c\k\a, where c is the phonon velocity 
given by 




- 2/3(n - 1)) - 2a-^7o- (22) 

In Fig. HI we plot the phonon velocity c as a function of the filling fraction n Q for 
a s /a = 0.01. We find that initially, for increasing uq, the phonon velocity increases. 
But for higher the phonon velocity starts decreasing due to the attractive effect of 
the occupation dependent tunneling term, until the phonon velocity becomes imaginary 
for a critical no- This results in a dynamical instability of the superfluid when we are 
within the limit < 1. This instability occurs due to the attractive effect of the 
occupation dependent tunneling, which can overcome the decreased repulsive on-site 
interaction depending on the number of particles per site no- To understand the effect 
of this instability, we first make a transition from the discrete Hubbard model to a 
continuous model applicable for fca < 1 with a continuous field <j>(r), 

#cont = - / d 3 r0*(r)V 2 0(O + ^ / V eS (r - OMOIW)! 2 - (23) 



2 „ 

Here, the distance is expressed with respect to the lattice constant a, and the 
effective interaction potential is given by V e s(r — r') = J r_1 [l + f3 — 0(no — 1) — 
2a^ 2 -7fc], where J 7-1 stands for the inverse Fourier transformation. Using a Gaussian 
ansatz along one direction, say x, and uniform in the other directions, <fr(r) = 
exp(—x 2 /2dl), the energy functional for the self-trapped state reads E so \ = 
\jd\ + -^%= (l + P- 2/?K - 1) - ct^(5 + exp(-2/d 2 s ))) ld s . When n exceeds a 
critical density, E so \ is minimized for a finite d s ^> 1. Thus, the homogeneous superfluid 
becomes dynamically unstable towards a state which is localized only in one direction, 
forming a 2D slab. 
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Figure 4. Phonon velocity c as a function of the superfiuid occupation number 
tiq. We find that after a critical occupation number, the phonon velocity becomes 
imaginary, denoting a dynamical instability. The fixed parameters are a s /a = 0.01 
and V /E r = 10. 

7. Conclusion and Outlook 

In this paper we have predicted various effects resulting from interaction-induced band 
mixing in systems of ultracold bosonic atoms in optical lattice potentials. We have 
derived the modified bosonic Hubbard model ((2]) having occupation-number-dependent 
parameters. This model comprises an effective interaction-induced broadening of the 
Wannier-like single-particle orbitals, and, thus, captures also the situation when the 
s-wave scattering length becomes comparable to the lattice spacing, a s /a — > 1. Using 
this model, we find that for scattering lengths a s ~ 0.15a and lattice depths Vq ~ 12Er, 
the n = 1 Mott-insulator state can become unstable towards a superfiuid which consists 
of bond-centered pair excitations. This scenario is novel considering the fact that the 
extended pairs emerge due to the occupation dependence of both the tunneling strength 
and the on-site interaction. For even higher interaction, the nature of the superfiuid 
pair excitations (destroying the insulator) changes. The pairs can now occupy both the 
bonds on the lattice (i.e., two neighboring sites) or its sites; in that way an exotic lattice 
geometry as shown in the center plot of Fig. [2] emerges. Increasing the interaction 
further, eventually the pairs live only on the sites of the lattice. In this regime of 
high interaction strength, the n = 2 Mott state gets completely destroyed by the 
site-centered pair fluctuations. We have also looked into the regime where interaction 
induced Wannier-broadening arises from large filling n > 1 at small scattering lengths, 
a s <C a. In this limit, we found that the superfiuid becomes dynamically unstable 
due to the attractive nature of the occupation-dependent tunneling. The system then 
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transforms from a uniform superfluid state to an asymmetric state which is localized in 
one direction and extended in the other two directions. 

In future studies, we would like to study the role of dissipation in these systems. 
Also, a more accurate determination of the number dependence of the Hubbard 
parameters J ni ,nj and U ni will be required for a quantitative description of the effects 
described here. Finally, it would also be worth studying in detail the role of a trapping 
potential, as it is present in experiments. 
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